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ABSTRACT 


It is not trivial to solve the equations that describe 
charged particle transport with the aid of computers, fo^ 
instabilities, inaccuracies, and subtle artifacts are well 
known afflictions of numerical analysis. Two specific points 
are: 


1. To avoid gross inaccuracies, pitch-angle scattering must 
be treated with great care. In particular, slightly inappro- 
priate numerical formulations give rise to mean free paths 
that are in error by large factors. 

2. A previously unrecognized artifact, "numerical disper- 
sion", is very similar to the physical phenomenon of dis- 
persion. To avoid misinterpretations arising from this 
similarity, the spatial increment of the finite-difference 
grid must be a small fraction of the mean free path. 

These points are illustrated by calculations based upon 
finite-difference approximations to the transport equation. 

1. Introduction . The diffusive idealization, which has been almost 
universally invoked in discussions of cosmic-ray transport is easy to 
treat analytically. However, many observed phenomena give clear evi- 
dence for the presence of non-dif fusive effects. One example is the so- 
called "scatter-free" propagation of kilovolt solar electrons which is 
inconsistent with diffusion, but which can readily be interpreted in 
terms of the coherent mode of propagation. Although the qualitative 
features of these effects have been outlined (Earl, 1974a, 1976), the 
theory is very complicated. Consequently, there is a need for reliable 
numerical computations which bypass these complexities and yield con- 
crete results that are well suited for comparison with observations. 
This paper explores such methods within the limited context of 
rectilinear propagation of cosmic-ray along a uniform guiding field on 
which are superimposed random fields. 


Under these circumstances, transport is described by 


a f ^ 3f _ 3 , af 

3s + ^3z 3p V 3p 


( 1 ) 


in which f is particle density in phase space, p is the pitch-angle co- 
sine, and z is distance parallel to the guiding field. The parameter 
s - Vt, where V is particle velocity, plays the role of a temporal vari- 
able. The coefficient of pitch-angle scattering is given by 


, _ (3/2 A) / * 2v , q-1 


v g (2-q)(4 ' -g ) < ‘ 1 " 1J ;I1J| ' ’ (2) 

where X is The mean free path, and q is an index that measures the 
anisotropy of scattering. 
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Diffusion refers to a configuration of slow temporal variations 
of the anisotropy and weak spatial inhomogenities of the isotropic den- 
sity F (Earl, 1974b), which act as the source of a diffusive anisotropy 
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where _ 


g{p} = [(4-q)X/3] Ulul 1 "* 1 (4) 0.8 


is a solution of 


±- * ii. = -y. 

9u 


q-1.5 


q = l.9 


2. Finite-Difference Approximations - 
to the Scattering Operator . In the g - 
discrete formulation, the continuous 

variables are replaced by a three- f , , , | , ■ i ■ i i i i 

dimensional grid whose spacings are -zp) inf) 

Az, Ay, and As = VAt, and the deriv- b IU 
atives appearing in equation (1) are Mmox 

replaced by their finite— difference , . 

analogs. A test of the relationships implied by diffusion can be made 
without extensive calculations by replacing the differential operator 
in equation (5) by its finite version and solving numerically for the 
diffusive anisotropy. This procedure should give an answer consistent 
with equation (4), and if it does not, the effective mean free path is 
given by 

X = (3/2) rj^ug dy. (6) 

Kota et al. (1982) gave an explicit finite-difference form for the scat- 
tering operator. Although this derivation seems eminently reasonable 
from the standpoint of both mathematics and physics, the resulting scat- 
tering operator does not pass this test near the threshold at q - 2 tor 
pure coherent propagation. 

Mean free paths derived from this operator are given in the 
figure above, where the ratio of the actual value of X obtained from 
equation (6) to the nominal value that appears in equation (2) is 
plotted against the number M = 2/ Ay of increments in pitch-angle. The 
correspondence is close for q = 1 and fairly close for q - 1 . 5 , but for 
q = 1 9 X is a factor of ~2.5 too small and shows little indication of 
converging to the nominal value as M increases. These inaccuracies are 
a consequence of the fact that weak scattering near y - 0 gives rise to 
large deviations from the analytic behavior. To avoid these deviations, 
I use a finite operator constructed so that equation (5) is satisfied 
exactly when g is described by the analytic expression. Far from y - 0, 
this operator behaves the same as Kota's operator, but it also g ves e 
correct diffusive aniso£tropy and makes the actual mean free path 
identical to the nominal one. 
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— — Numeri cal Dispersion . In a configuration of weak scattering, parti- 
cle bunches propagate coherently in a moving Gaussian profile given by 

F 0 - 1/ (4irD A t)^2 exp {-(z - V*t) 2 /4D*t> (7) 

where V* is a characteristic velocity, which corresponds to an angular 
distribution that is nearly isotropic in one hemisphere and zero in the 
other, and D* is the coefficient of dispersion, which describes the 
broadening of the Gaussian that arises from statistical fluctuations of 
individual velocities within this distribution. (Earl, 1974a). 

The finite— difference implementation leads to a very similar 
effect which can be described by equation (7) with D* replaced by 

D t = ^ ~ J As) (8) 

The total dispersion is a superposition of the physical effect described 
by D* and the numerical effect described by D f . To ensure that physics 
dominates, the condition D* » Df must be satisfied. In the examples 
presented below, where q = 1.8, this condition translates into a re- 
quirement that Az < A/40. This condition is similar in concept to the 
condition derived by Kota et al. for the validity of calculations of 
diffusion, which is Az < 3A/4. However, the required spatial resolution 
is much finer in the coherent regime that it is in the diffusion regime. 




— — A Specific Example. To illustrate the points made above, two calcu- 
lations were made which had identical values of A and Ay * 0.22 but 
which had values Az » A/21 (plotted as solid lines) and Az « A/210 
(plotted as dotted lines) that are, respectively, well below and well 
above the required spatial resolution. In the figure to the left 
snapshots of the density profiles at the instant s - 0.19 A, display the 
familiar coherent pulses followed by small but perceptible wakes. At 
this time, the dotted peak has moved a distance of 0.105 A to the right 
of the point indicated by a vertical dotted line at which a localized 
coherent angular distribution was injected at s - 0. This corresponds 
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to a coherent velocity ■ 0.55V that is in good agreement with the 
corresponding velocity of 0.53V tabulated by Bieber (1977). Similarly, 
the FWHM width of the pulse 0.135 X is in good agreement with the width 
0.125X expected from a combination of physical dispersion with 
D* = 0.0064XV and numerical dispersion with D^- 0.001 XV. 

The coherent velocity of the solid line peak is less than that of 
the dotted one and its dispersion, which is dominated by the numerical 
effect, is larger. To the right of point A, significant solid density 
appears in regions where the dotted density is zero. This discrepancy 
is relevant to descriptions of the early onset phase of solar particle 
events, for it means that the numerical dispersion can bring particles 
from the sun faster than the physics can, and this could lead to signi- 
ficant errors in the timing of flare events. 

The angular distributions in the figure that appears at the right 
above add weight to this point. They refer to point A where the two 
profiles have equal densities and, consequently, can be directly com- 
pared without normalization. The dotted distribution is much more 
collimated in the forward direction near y = 1 than the solid one is. 
It is clear that any attempts to interpret observed angular distribu- 
tions of the first few particles arriving in a flare event should take 
rigorous measures to exclude the effect of numerical dispersion. 

5. Conclusions . Before non-dif fusive transport can even be considered, 
the scattering operator must be carefully tailored to give valid results 
in the diffusion regime. Once this has been done, calculations that 
take proper account of numerical dispersion give results in good agree- 
ment with those expected from the theory of coherent propagation, but 
not with those calculated with inadequate spatial resolution. This sit- 
uation can lead to insidious errors, for the invalid results appear to 
be qualitatively reasonable, and can be detected only by quantitative 
tests . 

6. Acknowledgements . This research was supported by NASA under Grant 
NGR 21-002-066. 

References 

Bieber, J.W. 1977, Ph.D. thesis, University of Maryland. 

Earl, J. A. 1974a, Ap.J., 188, 379 (Paper II). 

Earl, J. A. 1974b, Ap.J., 193, 231 (Paper III). 

Earl, J. A. 1976, Ap.J., 205 , 900 (Paper IV). 

Kota, J., et al . 1982, Ap.J., 254 , 398. 



